Expression of microRNAs and their target genes in melanomas originating from gynecologic sites

Melanomas from gynecologic sites (MOGS) are rare and have poor survival. MicroRNAs (miRs) regulate gene expression and are dysregulated in cancer. We hypothesized that MOGS would display unique miR and mRNA expression profiles. The miR and mRNA expression profile in RNA from formalin fixed, paraffin embedded vaginal melanomas (relative to vaginal mucosa) and vulvar melanomas (relative to cutaneous melanoma) were measured with the Nanostring Human miRNA assay and Tumor Signaling mRNA assay. Differential patterns of expression were identified for 21 miRs in vaginal and 47 miRs in vulvar melanoma (fold change >2, p<0.01). In vaginal melanoma, miR-145-5p (tumor suppressor targeting TLR4, NRAS) was downregulated and miR-106a-5p, miR-17-5p, miR-20b-5p (members of miR-17-92 cluster) were upregulated. In vulvar melanoma, known tumor suppressors miR-200b-3p and miR-200a-3p were downregulated, and miR-20a-5p and miR-19b-3p, from the miR-17-92 cluster, were upregulated. Pathway analysis showed an enrichment of “proteoglycans in cancer”. Among differentially expressed mRNAs, topoisomerase IIα (TOP2A) was upregulated in both MOGS. Gene targets of dysregulated miRs were identified using publicly available databases and Pearson correlations. In vaginal melanoma, suppressor of cytokine signaling 3 (SOCS3) was downregulated, was a validated target of miR-19b-3p and miR-20a-5p and trended toward a significant inverse Pearson correlation with miR-19b-3p (p = 0.093). In vulvar melanoma, cyclin dependent kinase inhibitor 1A (CDKN1A) was downregulated, was the validated target of 22 upregulated miRs, and had a significant inverse Pearson correlation with miR-503-5p, miR-130a-3p, and miR-20a-5p (0.005 < p < 0.026). These findings support microRNAs as mediators of gene expression in MOGS.


Introduction
Melanoma is the fifth most common cancer in men and sixth most common cancer in women, accounting for 7 and 4% of all cancer diagnoses in each sex, respectively [1]. While melanoma most usually occurs in cutaneous sites, 6.8% of melanoma cases occur in non-cutaneous locations [2]. Mucosal surfaces are among the most commonly affected sites and account for up to 3.7% of all melanoma cases [3,4]. Women are at an almost two times greater risk of developing mucosal melanoma compared to men [5]. This gender discordance is influenced by the reported frequency of vulvar and vaginal melanoma, accounting for over 50% of mucosal melanomas in women [2].
Melanomas originating from gynecologic sites (MOGS) are a unique subset of melanoma tumors arising from regions of mucosa lining the female reproductive tract including the vagina and cervix, as well as the vulva. MOGS are rare, comprising only 1 to 3% of all melanoma cases diagnosed in women [6,7]. MOGS are reported to occur most frequently in the vulva and vagina, representing about 75% and 20% of mucosal melanoma cases, respectively [2]. Presentation of MOGS at the cervix is rare in comparison, encompassing only 3-9% of MOGS cases [8]. MOGS are associated with a poor clinical outcome and low survival rates, which may be due to a lack of accepted screening methods for early detection [8,9]. Patients with MOGS frequently present at an advanced disease stage due to their internal location. Pelvic and/or inguinal nodal involvement is an important prognostic factor and is reported in 25 to 50% of cases [10][11][12]. Surgical excision can be challenging given the proximity of affected tissues to anatomic structures of importance, such as the rectum and bladder [6,8,13]. Although novel, minimally invasive approaches are being utilized for gynecologic malignancies, such as neuropelveology, [14,15] MOGS frequently require vulvectomy, radical vaginectomy, inguinofemoral and/or pelvic lymphadenectomy and even total pelvic exenteration. However, the primary location of these tumors and aggressive nature of MOGS often require non-surgical management with systemic and/or local therapies. Unfortunately, there is a paucity of active agents, and novel targeted therapy is urgently needed to improve outcomes for this very rare and difficult to treat disease. Novel therapies must also minimize the psychological impact of treatment and preserve quality of life, which can be severely impaired during the treatment of gynecologic malignancies [16][17][18].
Previous reports have determined that MOGS have a lower frequency of oncogenic mutations in BRAF and NRAS than is observed in cutaneous melanoma not associated with chronic UV-damage [10]. Sequencing of frequently mutated oncogenes associated with melanoma demonstrates that these variants are infrequently observed in MOGS. In MOGS, BRAF and NRAS mutations occur at an incidence of 2.6% and 5.3%, respectively, compared to 59-65% and 20% in cutaneous melanoma [7,10,19]. c-KIT mutations are observed with slightly higher frequency in mucosal melanomas with an estimated incidence of 22.2 to 25%. The use of KIT inhibitors such as imatinib has resulted in only a 16% response rate among these patients [7,10]. Thus, there is a need for characterization of MOGS beyond mutational status of known oncogenes found in cutaneous melanoma to determine the unique molecular features that may contribute to MOGS progression.
The rarity of MOGS has hindered development of a consensus regarding the standard management of these tumors [6]. Some phase III trials using ipilimumab or BRAF V600E targeting agents have excluded patients with primary mucosal or vulvovaginal melanomas [10]. One study that compared survival of patients with mucosal melanoma receiving immunotherapy to cutaneous melanoma receiving immunotherapy found that patients with vulvovaginal melanoma had a median overall survival of 8.6 months, while the median overall survival in patients with cutaneous melanoma was 14.5 months [20]. Additionally, some trials that did not specifically exclude these patients did not collect details pertaining to the primary tumor site [10]. Therefore, the potential response of MOGS to current cutaneous melanoma therapies, as well as the correlation of gene expression with response to therapy, has not been thoroughly assessed [21]. microRNAs (miRs) are non-coding RNA sequences 20 to 23 nucleotides in length that function to regulate expression of gene targets by binding with mRNA transcripts, resulting in mRNA degradation and/or translational inhibition [22]. miR expression profiling by our group and others has demonstrated that specific miRs are overexpressed in cutaneous melanocytic lesions and metastatic melanoma and can contribute to disease progression. miR profiling could be of prognostic and diagnostic value, including distinguishing MOGS from other gynecologic lesions [23][24][25][26]. Additionally, therapeutic approaches to modulate miRs in cancer are currently being developed [27,28]. However, a review of the literature reveals no reports evaluating miR expression specifically in MOGS [23,24,[29][30][31][32][33][34][35]. Therefore, the role of miRs in the development and progression of MOGS is unknown.
We aimed to determine the miR expression profile of MOGS and its potential effects on tumor signal transduction and gene expression. In the present study, miR and mRNA expression profiles were analyzed in melanomas originating from the vagina relative to patient matched, normal adjacent vaginal mucosal tissue using the Nanostring platform [36]. Relative miR and mRNA expression patterns in vulvar melanoma tissue were determined relative to primary cutaneous, nongynecologic melanoma as a means of evaluating the importance of the comparator tissue. The correlation between miR expression and mRNA expression of predicted and known target genes and the potential functional impact of their interaction in MOGS were then examined.

Sample selection and tissue collection
6 samples of vaginal melanoma with paired vaginal mucosa and 22 samples of vulvar melanoma tissue with 9 samples of primary cutaneous melanoma were included in the analysis. De-identified samples of formalin fixed, paraffin embedded (FFPE) vaginal melanoma and adjacent normal vaginal mucosal tissue, as well as vulvar melanoma tissue were selected from the tissue archive at the University of Virginia for use under an approved IRB protocol no. 2007C0054. FFPE samples of vulvar melanoma and primary cutaneous melanoma were also obtained from the pathology core facility tissue archive at the Ohio State University under IRB protocol no. 2007C0015. Melanoma tissue storage was registered with ClinicalTrials.gov (Clin-icalTrials.gov Identifier: NCT04567706).
Eligible patients were those that had histologically confirmed vaginal or vulvar melanoma treated with surgery prior to receiving systemic therapy. Given that the goal of the study was to understand overall MOGS biology, no other criteria were applied. For vaginal melanoma samples, 2 mm diameter cylindrical punch samples of the FFPE tumor tissue and adjacent vaginal mucosa were isolated from the paraffin blocks using sterile technique. Regions of interest selected for tissue collection were identified based on review of hematoxylin and eosin-stained tissue sections by a dermatopathologist (CC). For vulvar melanoma and primary cutaneous melanoma samples, four, 20 μm thick scrolls were collected from the tissue block for each sample. Biopsy samples were identified for inclusion based on review of hematoxylin and eosinstained sections by a dermatopathologist (CC).

RNA isolation
Two methods of RNA isolation were employed to optimize use of available and appropriate control tissues for comparison of miR expression in each tumor site. For vaginal melanoma and adjacent mucosal tissues, RNA was isolated from punch samples of formalin fixed, paraffin embedded (FFPE) melanoma or normal adjacent vaginal mucosal tissue using a modified Qiagen miRNeasy FFPE kit protocol (Qiagen, Hilden, Germany). Briefly, paraffin was washed from samples of FFPE tissue with xylene and 100% ethanol prior to Proteinase K digestion overnight at 50˚C with constant agitation. Samples were then incubated at 80˚C, then centrifuged at maximum speed. DNA digestion was then performed on the supernatant. The resulting solution was then added to an RNeasy MinElute spin column, washed and dried by centrifugation prior to RNA elution in nuclease-free water. Isolated RNA was stored at -80˚C and assessed for concentration and purity by Nanodrop and Qubit analysis.
For vulvar melanoma and primary cutaneous melanoma samples, RNA was isolated from FFPE tissue scrolls using the Invitrogen RecoverAll Total Nucleic Acid Isolation Kit according to manufacturer's instructions (Invitrogen, Carlsbad, CA). Briefly, four 20 μm-thick tissue scrolls were washed with xylene and 100% ethanol prior to protease and DNAse digestion and sample elution. The isolated RNA was then further concentrated using the Norgen RNA Clean up and Concentration kit according to manufacturer's instructions (Norgen, Thorold, ON, Canada). All RNA was stored at -80˚C and assessed for concentration and purity by Nanodrop and Qubit analysis.

NanoString microRNA expression assay
Isolated and purified RNA (100 ng) was loaded onto a NanoString nCounter (NanoString Technologies, Seattle, WA) cartridge and quantification of miR expression was carried out as previously described [37]. The fluorescence of each hybridized miR-probe complex was analyzed by an nCounter Digital analyzer (NanoString Technologies, Seattle, WA) via high-density scan with 600 fields of view. The investigation included 5 positive, 5 negative, and 5 housekeeping genes provided by the manufacturer.

Nanostring microRNA assay analysis
All vaginal melanoma and vaginal mucosa samples were included on a single cartridge. For vulvar melanoma, both vulvar melanoma and cutaneous melanoma samples were included on each cartridge and distributed across three cartridges. Negative control miR targets included in the panel were used to assess background hybridization and filter out target miRs with low expression. miRs with mean expression levels lower than the highest detected negative control sample were removed from the analysis. The upper quartile normalization method was then used to normalize across biological samples. Differential expression of miRs between vaginal melanoma and normal vaginal mucosa or vulvar melanoma and primary cutaneous melanoma was detected using linear models and a moderated t-test while considering variation from cartridge effect [38]. Significance was adjusted using a Bonferroni procedure to control the mean number of false positives [39]. Statistical software SAS 9.4 and R 3.6 was used for analysis, with an α of 0.01 and fold change of at least 2 used to identify differentially expressed miRs.

Pathway enrichment analysis for differentially expressed miRs
Pathway enrichment analysis was performed to assess the impact of miRs on canonical signaling pathways [40]. The gene targets of miRs with a significantly decreased expression were identified using the DIANATools miRPath database V.3 (DIANA-mirPath, RRID: SCR_017354) (http://snf-515788.vm.okeanos.grnet.gr/). These genes were then used to identify which Kyoto Encyclopedia of Genes and Genomes (KEGG) pathways may be impacted by downregulated miRs [40,41]. Pathways significantly enriched by miRs with decreased expression in vaginal or vulvar tumors were identified based on experimentally validated interactions listed in TarBase, using Fisher's Exact Test (hypergeometric distribution, p-value threshold = 0.05). The Benjamini-Hochberg procedure was used to control error rate (false discovery rate) due to multiple comparisons. This procedure for KEGG pathway enrichment was repeated for upregulated miRs in vaginal and vulvar melanoma [42].

Nanostring mRNA expression assay
Expression of 780 mRNA transcripts and 20 housekeeping genes was then assessed in samples of RNA isolated from vaginal melanomas and normal vaginal mucosa, or vulvar melanoma and primary cutaneous melanoma using the Nanostring nCounter Tumor Signaling 360 assay (Nanostring Technologies, Seattle, WA). RNA concentration and purity was confirmed using both Nanodrop spectrophotometry and Qubit High Sensitive RNA Assay (Thermo Fisher Scientific). Normalized quantitation of RNA for loading into the nCounter1 system was determined using the percentage of RNA fragments that are larger than 200 nucleotides in size (DV200), using a 2200 Tapestation (Agilent) with an RNA Screentape kit. For DV200 > 65% 100 ng of RNA was used; for DV200 between 25% and 65%, 200 ng of RNA was used; for DV200 < 25%, 300 ng of RNA was used. Quantification of mRNA expression was performed according to the manufacturer's instructions.
Two standard control RNA samples supplied by the manufacturer were included in the assay to permit normalization of mRNA expression in samples across multiple cartridges. Six positive and eight negative control targets were included on each cartridge to assess background hybridization and remove target mRNAs with low expression. The 780 mRNA targets included in this assay are as defined by the manufacturer. All mRNAs that demonstrated mean expression below that of the highest detected negative control samples were removed from the analysis. The upper quartile normalization method was used to normalize across biological samples. Differentially expressed mRNAs were identified as those demonstrating a minimum of 2-fold change in expression with an adjusted p-value of 0.01, determined using a t-test followed by a Bonferroni correction for multiple comparisons.

Pathway enrichment analysis and functional annotation for differentially expressed mRNA transcripts in gynecologic melanoma tissue
Ingenuity Pathway Analysis, RRID:SCR_008653 (IPA, Qiagen) was employed for identification of the canonical signaling pathways that were significantly enriched by differentially expressed mRNA transcripts in vaginal and vulvar melanoma tissue. This platform was also used to identify the functional annotations for each differentially expressed gene. Both pathway enrichment and functional annotation enrichment were performed using the Core Analysis function in IPA (p-value cutoff = 0.05). The reference gene set used for both canonical pathway enrichment and functional annotation was comprised of all genes assessed for mRNA expression in the Nanostring nCounter Tumor Signaling 360 panel to minimize enrichment of false positive annotations and pathways.

Identification of differentially expressed mRNAs with predicted regulation by multiple differentially expressed microRNAs
Differentially expressed miRs were determined based on mean normalized expression across all vaginal or vulvar melanoma samples relative to normal vaginal mucosa or primary cutaneous melanoma samples. Experimentally validated mRNA targets for these differentially expressed miRs were identified using the miRNet database (https://www.mirnet.ca/) and this list was further restricted to include only those mRNAs that were: 1) confirmed via the present Nanostring analysis; 2) inversely correlated in the analysis; and 3) were the target of at least two differentially expressed miRs [43]. A miR-gene target interaction network was then generated in miRNet for each set of downregulated miRs and their corresponding upregulated target genes in vaginal and vulvar melanoma. This procedure was repeated to generate separate miRgene target interaction networks for significantly upregulated miRs and downregulated target genes [42].
Pairwise correlation analysis of microRNA and mRNA expression within melanomas originating from gynecologic sites An alternative bioinformatics approach was also employed in order to explore potential one-toone miR-mRNA interactions within and across individual tumor samples. Pearson correlation analysis was performed between each differentially expressed miR and each differentially expressed mRNA transcript in vaginal or vulvar melanoma. Expression of each possible miR-mRNA pair was assessed within all individual tissue samples for evidence of a significant inverse correlation in expression within the group that would support a direct inhibitory effect of the miR on expression of the target gene mRNA transcript. A Pearson correlation coefficient ("r", ranging from -1 to 1) was determined for each pair, with those less than 0 considered consistent with an inverse correlation in expression. From this list of miR-mRNA pairs with a significant inverse correlation, the microRNA target filter function in the Ingenuity Pathway Analysis, RRID:SCR_008653 (IPA) software (Qiagen) was used to select for predicted (based on TargetScan database) or experimentally validated (based on miRTarBase database) target genes for each differentially expressed miR. The strength of the correlation for each miR-mRNA pair was determined by the Pearson correlation coefficient, as described above. Significance of each Pearson correlation was determined using a t-test for linear regression, based on an α of 0.05 [44].
In vulvar melanoma, miRs with decreased expression relative to primary cutaneous melanoma resulted in significant enrichment of 11 pathways, among which "microRNAs in cancer" was most significantly enriched by 2 miRs (miR-200a-3p and miR-200b-3p, p = 0.011 and p = 1.83E-09, respectively) based on interaction with 43 validated target genes in the pathway (Fig 2B). For miRs with increased expression in vulvar melanoma relative to primary cutaneous melanoma, there was an enrichment of 35 pathways, among which "proteoglycans in cancer" was again the most significantly enriched by 26 miRs (6.0E-12 < p < 0.041) based on the interaction with 153 validated target genes known to contribute to the pathway (Fig 2C).

Differential expression of mRNA transcripts in melanomas originating from gynecologic sites
In order to further characterize the impact of microRNA expression patterns on gynecologic melanoma biology, mRNA expression was evaluated. Differentially expressed mRNA transcripts were identified in vaginal melanoma tissue relative to normal vaginal mucosa, as well as in vulvar melanoma when compared to primary cutaneous melanoma tissue using the Nanostring nCounter Tumor Signaling 360 assay for mRNA expression. In vaginal melanoma 10 mRNA transcripts had significantly increased expression including TPX2, FOXM1,   Fig 3A). All differentially expressed mRNA transcripts in vaginal melanoma relative to normal vaginal mucosa are listed in Table 5 and functional annotations of these genes are listed in Table 6. In vulvar melanoma, 89 genes exhibited differential patterns of expression on the mRNA level relative to primary cutaneous melanoma. 43 mRNA transcripts were expressed at significantly higher levels in vulvar melanoma tissue with DLL3, NUF2, AUKRA, RFC3, and TOP2A being among the most highly upregulated transcripts. 46 mRNA transcripts were expressed at significantly lower levels in vulvar melanoma tissue with KRT17, CALML3, FGFR2, TPSAB1/

PLOS ONE
microRNA expression in melanomas originated from gynecologic sites B2, and SERPINB5 being among the most downregulated transcripts (fold change > 2, p < 0.01 for all, Fig 4A). All differentially expressed mRNA transcripts in vulvar melanoma relative to primary cutaneous melanoma are listed in Table 7 and functional annotations of these genes are listed in Table 8.

Pathway enrichment analysis for differentially expressed mRNA transcripts in melanomas originating from gynecologic sites
Pathway enrichment analysis was then completed using the Core Analysis function in IPA for all differentially expressed mRNA transcripts in vaginal tissue relative to vaginal mucosa and

PLOS ONE
microRNA expression in melanomas originated from gynecologic sites vulvar melanoma relative to primary cutaneous melanoma. This analysis permits identification of cellular signaling pathways and associated functions impacted by the differential pattern of gene expression between groups. In vaginal melanoma 8 pathways were significantly enriched by the previously identified differential expression of 17 mRNA transcripts (10 up, 7 down) relative to normal vaginal mucosa (p < 0.05, Fig 3B). Of these, the most significantly enriched pathway was "RhoA signaling" (p = 4.59E-03) which incorporates 2 of the mRNA transcripts identified in vaginal melanoma with significantly decreased expression relative to normal  Fig 4B). Of these, the most significantly enriched pathway was "Kinetochore metaphase signaling pathway" (p = 5.38E-05) which incorporates 9 genes with significantly increased mRNA expression in vulvar melanoma (AURKB, BIRC5, BUB1, BUB1B, CCNB1, CDK1, KIF2C, NUF2, and PLK1).

Differentially expressed mRNA transcripts targeted by two or more differentially expressed miRs in melanomas originating from gynecologic sites
In order to identify mRNA transcripts regulated by multiple miRs, miR-mRNA interaction networks were constructed using validated miR-mRNA interactions listed in TarBase. Given the inhibitory nature of miR activity, reduced miR levels would be expected to result in increased mRNA expression and vice versa. In vaginal melanoma, MCM4 was the only gene with increased mRNA expression that was identified as a target of two more downregulated miRs: miR-99a-5p and let-7c-5p ( Fig 5A). Conversely, FOS and SOCS3 were both downregulated in vaginal melanoma and were identified as validated gene targets of two or more upregulated miRs (miR-17-5p, miR-181a-5p, and miR-20a-5p for FOS, and miR-19b-3p and miR-20a-5p for SOCS3, Fig 5B).
In vulvar melanoma, mRNA transcripts for BUB1, NME1, RRM2, CCNE1, and KPNA2 were upregulated and identified as target genes for two or more miRs with significantly decreased mean expression, including miR-494-3p, miR-200b-3p, and miR-200a-3p (Fig 6A). Among the 46 genes with decreased mRNA expression in vulvar melanoma, 25 were identified as targets of two or more upregulated miRs (Fig 6B). The three genes associated with the greatest number of dysregulated miRs were CDKN1A, AXL, and EGFR, which are validated targets for 22, 14, and 7 upregulated miRs, respectively.

microRNA-mRNA pairs with an inverse Pearson correlation and predicted or validated interaction in expression in melanomas originating from gynecologic sites
Individual miR-mRNA pairs were then assessed by Pearson correlation analysis for evidence of an inverse correlation in expression within individual tissue samples and maintained across all samples. For vaginal melanoma, a total of 17 unique miR-mRNA pairs were identified with a negative correlation (r < 0) and a p-value less than 0.05 (Fig 5C). However, none of the 17 miR-mRNA pairs exhibited predicted or validated interactions. The analysis was repeated to Table 8. Significantly enriched cellular and molecular functions for differentially expressed genes in vulvar melanoma.

Discussion
The rarity of MOGS has resulted in poor molecular characterization, particularly relative to cutaneous melanoma. In the current study, microRNA and mRNA expression were characterized in vaginal and vulvar melanoma using high throughput assays for assessment of hundreds of miR and mRNA targets per sample. The potential impact of miR and mRNA expression in MOGS was examined using two unique approaches. These included identification of genes with potential regulation by multiple predicted miRs and also an assessment of 1-to-1 Pearson correlations between individual miRs and target genes across individual samples. Together, these two approaches permitted the identification of the miR-mRNA pairs that have the greatest potential to interact in a functional manner in the setting of MOGS. Among the 5 downregulated miRs in vaginal melanoma, miR-145-5p has been shown to function as a tumor suppressor miR in cutaneous melanoma via inhibition of TLR4 and NRAS expression and subsequent regulation of the downstream NF-κB and MAPK signaling pathways that contribute to tumor proliferation, migration, and invasion [45,46]. There were 14 significantly upregulated miRs in vaginal melanoma and miR-106a-5p, miR-17-5p, and miR-20b-5p displayed the greatest fold change. Elevated expression of miR-17-5p and miR-106a-5p has previously been detected in the plasma of stage III-IV melanoma patients by our group [47]. miR-17-5p, a member of the miR-17-92 microRNA cluster, has been described as an oncomiR in melanoma based on its ability to enhance cell proliferation in patient-derived primary melanoma cultures [48]. In a similar fashion, miR-106a-5p and miR-20b-5p are both members of the miR-106a-363 cluster, which has been associated with cell cycle progression in melanoma [49]. miR-579-3p displayed a smaller increase in vaginal melanoma (log 2 (fold change) = 2.7368). Low expression of miR-579-3p has been previously reported to be a negative prognostic marker in cutaneous melanoma, where decreased expression has been hypothesized to promote resistance to BRAF and MEK inhibitors [50]. Here we report the opposite finding (increased expression in vaginal melanoma), which may reflect the intrinsic differences between cutaneous melanoma and MOGS.
Of the 3 miRs with significantly decreased expression in vulvar melanoma, two are of interest, namely miR-200b-3p and miR-200a-3p. miR-200b and miR-200a have been shown to function as tumor suppressor miRs in melanoma [51,52]. Thus, it is not surprising to learn that their expression is directly inhibited in melanoma cell lines by an overexpressed long noncoding RNA lncRNA-HEIH that binds to their promoter [52]. This inhibition of miR-200b and miR-200a has been shown to result in increased proliferation, migration, and invasion of melanoma cells. Among the 44 miRs with increased expression in vulvar melanoma, miR-20a-5p+20b-5p and miR-19b-3p had the greatest increase in expression. miR-20a-5p and miR-19b-3p are members of the miR-17-92 cluster mentioned above, and their upregulation suggests the involvement of this miR group in pathogenesis of both vulvar and vaginal melanoma. Indeed, an increase in expression of these two miRs contributed to the significant enrichment of the "proteoglycans in cancer" pathway, which was the most significantly enriched pathway among miRs with increased expression in both vaginal and vulvar melanoma. Proteoglycans provide a contact link between the cell membrane and the surrounding extracellular matrix and play a central role in role in cancer cell adhesion and migration [53]. Descriptions of miR-17-92 cluster members as contributors to melanoma tumorigenesis or progression have been limited in cutaneous melanoma, but there are reports of its role in non-cutaneous cancers [54].
In order to characterize the impact of microRNA expression patterns in gynecologic melanoma, mRNA expression was also evaluated. Topoisomerase IIα (TOP2A) had significantly increased mRNA expression in both vaginal and vulvar melanoma relative to normal vaginal mucosa and primary cutaneous melanoma, respectively. TOP2A is critical for accurate chromosome segregation during mitosis and is commonly altered at the genome level in a variety of human cancers [55]. Significantly increased expression of TOP2A was also previously identified in vulvovaginal melanoma by immunohistochemistry in a study of 51 cases relative to 2253 non-gynecologic melanomas [56]. Several other genes with increased expression in vulvar melanoma have not been previously described as being upregulated in this setting. Examples are DLL3, NUF2, AURKA, and RFC3. NUF2 and AURKA are both involved in cell division and both contributed to the enrichment of the top four functional categories in vulvar melanoma ("cell cycle", "cellular assembly and organization", "DNA replication, recombination, and repair", and "cell morphology") ( Table 6).
In order to determine how differentially expressed (DE) miRs may function together, experimentally validated targets of multiple miRs were evaluated using miRNet. Through this approach, MCM4 was identified as an upregulated (log2(fold change = 2.9789), validated target gene of two significantly downregulated miRs in vaginal melanoma (let-7c-5p and miR-99a-5p) (Fig 5A and 5B). The MCM4 gene (minichromosome maintenance protein 4) codes for a protein with helicase activity that, when assembled into an MCM protein complex, functions in DNA replication [57]. MCM4 expression also contributed to the significant enrichment of two of the top three functional annotations in vaginal melanoma based on gene expression: "cellular growth and proliferation", and "DNA replication, recombination, and repair". Increased MCM4 expression has a significant association with reduced overall survival in both cutaneous melanoma and non-melanoma tumors arising from gynecologic sites [58,59]. The prognostic significance of increased expression of MCM4 is worthy of further study.
FOS and SOCS3 had significantly decreased expression in vaginal melanoma tissue (fold change = -3.121 and -2.7946, p = 0.00279 and 0.00978) and are validated target genes of multiple miRs that were significantly upregulated in vaginal melanoma. SOCS3 also had an inverse trending Pearson correlation across individual samples with miR-19b-3p (p = 0.093, r = -0.73). SOCS3 (suppressor of cytokine signaling 3) is suspected to have an anti-oncogenic role via its ability to suppress STAT3 activation in the JAK-STAT pathway that mediates signal transduction for the interferons [60].
In vulvar melanoma, this same approach identified 5 upregulated genes (BUB1, NME1, RRM2, CCNE1, and KPNA2) that are validated targets of downregulated miRs, with each target gene regulated by two of the three downregulated miRs in vulvar tumors. Among these genes, NME1 is a validated target of miR-200a-3p and miR-494-3p (both significantly downregulated) (Fig 7B). NME1 functions as a driver of cell proliferation and expression of stem-like markers in melanoma cell lines when expression is high [61]. CCNE1, a stem-like marker, was also elevated in vulvar melanoma relative to primary cutaneous melanoma. Therefore, loss of miR-200a-3p and miR-494-3p together with increased expression of their target gene NME1 may play a role in MOGS stemness [61].
Conversely, 25 downregulated genes were identified as experimentally validated targets of two or more miRs with increased expression in vulvar melanoma. Among them, CDKN1A (cyclin dependent kinase inhibitor 1A) was significantly downregulated (log 2 (fold change) = -1.173) and found to be a validated target gene for 22 miRs with increased expression in vulvar melanoma. CDKN1A also specifically demonstrated a significant inverse Pearson correlation in expression across individual samples with three miRs predicted or validated to regulate its expression: miR-503-5p, miR-20a-5p, and miR-130a-3p. CDKN1A encodes for p21, a cyclin dependent kinase regulator for the G1/S checkpoint of the cell cycle. p21 can function as either a tumor suppressor or an oncogene depending on subcellular localization of the protein [62]. Of particular importance in melanoma is the role of p21 (CDKN1A) as a principal effector of p53-mediated cell senescence.
The current study is the first to examine the unique miR expression profile MOGS. The use of two bioinformatics approaches for data analysis permitted the identification of the miR-mRNA pairs that have the greatest potential to interact in a functional manner in the setting of MOGS. A few limitations are noted. The small sample sizes used in this study are a function of the rarity of these tumors. Inter-chip variability for assessment of microRNA and mRNA expression in vulvar melanoma is also considered a limitation of the Nanostring platform but was limited through the inclusion of replicate samples on each chip. Additionally, the tissues used for comparison to tumor tissue expression differed between groups. In particular, vaginal melanoma was compared to adjacent healthy vaginal mucosa, while vulvar melanoma was compared to primary cutaneous, non-gynecologic melanoma as a means of comparing the importance of the comparator tissue. While there was some overlap in the microRNAs identified when reference groups were switched, additional work is needed to identify the best control tissue to characterize these rare types of melanoma. Additional work is also needed to determine if the miR profiles identified here correlate with survival or are predictive of therapeutic response. Due to the use of bulk tissue RNA, spatial heterogeneity of miR and mRNA expression cannot be determined from the current study. Lastly, the current study is limited to the expression of miRs and their target genes at the mRNA level only, and therefore confirmation of changes at the protein level will be pursued in future studies.
The current study has demonstrated unique patterns of miR and mRNA expression associated with both vaginal and vulvar melanoma. The identification of both novel and welldescribed genes and their oncogenic pathways provides targets of potential therapeutic utility for future investigation in the treatment of MOGS. miR profiles could thus be employed in the characterization of gynecologic pigmented lesions. They could also be used in the setting of socalled liquid biopsies where plasma levels of a miR could help to make the diagnosis of MOGS or follow the response to systemic therapy. The identification of miR-mRNA pathways common to MOGS in multiple patients could provide a rationale for future therapeutic approaches with targeted agents. Also, miR antagonists are available that could be delivered either systemically or intra-tumorally in an experimental setting. However, additional work is needed to further investigate the specific miRs and miR-mRNA interactions identified here to evaluate their potential therapeutic utility in the setting of MOGS.